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ABSTRACT 

We investigate the far-infrared-radio correlation (FRC) of stellar-mass-selected galax- 
ies in the Extended Chandra Deep Field South using far-infrared imaging from Spitzer 
and radio imaging from the Very Large Array and Giant Metre- Wave Radio Telescope. 
We stack in redshift bins to probe galaxies below the noise and confusion limits. Radio 
fluxes are -fC-corrected using observed flux ratios, leading to tentative evidence for an 
evolution in spectral index. We compare spectral energy distribution (SED) templates 
of local galaxies for ii'-correcting FIR fluxes, and show that the data are best fit by 
a quiescent spiral template (M51) rather than a warm starburst (M82) or ULIRG 
(Arp220), implying a predominance of cold dust in massive galaxies at high redshift. 
In contrast we measure total infrared luminosities that are consistent with high star- 
formation rates. We observe that the FRC index (g) does not evolve significantly over 
z = — 2 when computed from ii'-corrcctcd 24 or 160-/im photometry, but that using 
70- /zm fluxes leads to an apparent decline in q beyond z ~ 1. This suggests some 
change in the SED at high redshift, either a steepening of the spectrum at rest-frame 
~ 25 — 35/j,m or a deficiency at '^ 70^m leading to a drop in the total infrared/radio 
ratios. We compare our results to other work in the literature and find synergies with 
recent findings on the high- redshift FRC, high specific star-formation rates of massive 
galaxies and the cold dust temperatures in these galaxies. 

Key words: galaxies: high-redshift - galaxies: evolution - galaxies: ISM - infrared: 
galaxies - radio continuum: galaxies. 



1 INTRODUCTION 

One of tlie most exciting research areas in observational 
astronomy at this time is the field of far-infrared astron- 
omy. Our understanding of extragalactic sources in the 
far-infrared (FIR) and sub-millimetre (sub-mm) to mil- 
Umetre regimes has improved exceptionally over the past 
decade, thanks to such instruments as E SA's Infrared 
Spac e Observatory {ISO, launched 1995; iKessler et al..l 
19961') and NASA's S pitzer Space Telescope (launched 2003; 
Werner et al.J [200J), alongside ground-based instruments 
such as th e Sub-millimetre Com mon-User Bolometer Array 
(SCUBA; iHoUand et Z1 1 19991 ) and t he Max Planck Mi l- 
Umeter Bolometer Array (MAMBO: iKrevsa et al.J 1 19981 ). 
both commissioned in the late 1990's. One of the most im- 
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portant FIR instruments to dat e is the Mult i band Imaging 
Photometer for Spitzer (MIPS: Ipleke et al..|[20o3 '). which, 
with Spitzer's 0.85m mirror provides imaging with diffrac- 
tion limited resolutions of 6, 18 and 40 arcsec in three broad 
bands centred at 24, 70 and 160^m respectively. 

Recent advances in the resolution of the cosmic infrared 
bac kground (CIB) by stacking into images from Spitzer (see 
e.g. iDole et al. .112000 : iDve et al.|[2007l : ICharv fc Pop3l2010D 
and BLAS T (The Balloon-born e Large Aperture Sub-mm 
Telescope; iMarsden et al.j|2009l ') have enabled an improved 
understanding of the history of star formation and galaxy 
formation and evolution, and this will be further improved 
by ongoing work (such as lBerta et al..ll2010l) with E SA's new 
Herschel Space Observatory ( Pilbratt et al.j|2010l ). 

Most of the stellar mass in the local universe is con- 
centrated in the most massive galaxies (A/* > lO^^M©; 
iKauffmann et al..l 120031 ) and observations show that these 
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have been in place sinc e 



(IDickinson et al.1 



Bundv. Ellis fc Conselicd 2003: Perez- Gonzalez et al. 



2002; 



Taylor et al.j|2009l : ICoUins et al.j|2009l ). The significant 



2008 



crease in density of luminous (Lg_iooonm > IO^^Lq) and 
ultra- luminous (Z/g-iooo/jm ^ lO^^L©) infrared galaxies 
(LIR GS and ULIRGs) fro m the local universe to 2 ~ 2 — 3 
(e.g. iDaddi et al.1 l2005al lbl: ICaputi et all I2OO6I ) is thought 
to reveal the formation stages of these latter-day giants, 
which apparently formed in a remarkably short time between 
1 S; -2 ^ 3, in an antithe sis to the paradigm of hierarchical 
structure formation (e.g. Ide Lucia et al.ll2006l ). 

In the FIR-sub-mm regime the dominant source of the 
extragalactic background light (after the cosmic microwave 
background) is thermal continuum emission from interstel- 
lar dust, which is mainly composed of polycyclic aromatic 
hydrocarbons (PAHs), gr aphites and silicates , typically less 
than ~ 0.25/xm in size (jDraine fc Lil 120071 : iDraine et al..l 
[20071). The FIR emission in star- forming galaxies is thought 
to arise both from cold dust in the large-scale 'cirrus' com- 
ponent of the interstellar medium (ISM), an d from warmer 
dust i n and around star-forming regions fe.g. Ide Jong et al.l 
1 19841 : iHeloulllOsd ). If there is a sufHcient level of dust- 
enshrouded star formation then the FIR emission is dom- 
inated by this 'warm' dust component , which has character- 
istic temperatures of around 30-50K (IDunne fc Bales I2OOII : 



ISaiina et ai]l2006l : IPve et al.ll2007l : IPascale et al.. 1120091 ) .The 
dust is heated to these temperatures by the ultraviolet (UV) 
radiation field from hot O and B type stars, which are 
present only in regions of ongoing or recent star forma- 
tion ('recent' m eaning within the lifetime of these short-lived 
stars, < 10^ vr: lKennicuttlll998l '). For this reason, the total 
IR luminosity (Ltir = is-iooo/im) can be used as a tracer 
of star-formation rates (SFRs) in galaxies, often using one 
or a number of FIR fluxes (such as t he Spitzer MIPS bands) 
or sub- mm fluxes to estimate Ltir (|Kennicuttll 19981 ) . 

One issue with using the FIR as a SFR tracer is the 
contribution from cold dust in the ISM, which is heated by 
older stars in the disk of the galaxy, and is th erefore un- 
related to star formation (jCalzetti et al..ll20ld V Contami- 
nation of samples by galaxies hosting active galactic nuclei 
(AGN) is another problem, as AGN also emit UV radia- 
tion which can heat dust in the torus. AGN-heated dust is 
generally hotter than dust heated in star-forming regions, 
so the thermal spectrum peaks at a shorter wavelength, 
and mid-infrared (MIR) fluxes (including Spiteer's 24-/im 
band, as well as shorter wavelengths) would be boosted. 
Mid-infrared fluxes are also affected more uncertainly by the 
emission features of PAH molecules, which are ubiquitous in 
star-forming galaxies (e . g. iLeger fc Pugedll984l:lRoche et al 



I991I : iLutz et al.1 1 19981 : lAUamandola. Hudgins fc Sandford 



19991 ). and the 10- /.im silicate absorption trough, so the 



use of 24-^m fluxes as SFR indicators at high redshifts 
is subject to s ome c ontention (se e e.g. iDale et al. 20051: 



Calzetti et al..l l2007l: iDaddi et al 



2OO7I : lYoung. Bendo fc LucerdbOOgJ : bieke et al.ll2009h 



20071: iPapovich et al.] 



Another part of a galaxy's SED that can be used as 
a SFR indicator is the radio luminosity. Non-thermal ra- 
dio continuum emission from star-forming galaxies origi- 
nates from type II supernova remnants (SNRs), the end- 
points of the same massive short-lived stars that heat 
the dust via their UV radiation. This connection with 
dust heating is important because it leads to the well- 



known (but not fully underst o od) FIR-Radio correla- 
tion (FRC: Ivan der KruitI Il973l: iRickard fc HarvevI [l98i : 



iHelou. Soifer fc Rowan-Robinsonill985l : ICondonlll992l . etc). 
The FRC is linear, remarkably tight and holds for a wide 
range of galax y types over at least five orde rs of magnitude 
in luminosity (|Yun. Reddv fc Condonll200lh . It can be ex- 
plained in terms of ongoing star formation producing hot 
massive (M > 8M0) stars: while the FIR flux is emitted 
from dust heated by these stars, the radio emission arises 
from synchrotron radiation by cosmic ray (CR) electrons ac- 
celerated in the SNRs of the dying stars. The non-thermal 
radio emission is smeared out through the galaxy as the 
relativistic CR electrons travel through the galaxy over life- 
times of r ~ 10* years, during which they emit synchrotron 
radiation via int eractions with the galactic magnetic field 
(a s described by ICondonlll992l . and shown observationally 
bv lMurphv et al.. 1120061 ). A shallower thermal component is 
also present in the radio spectrum due to bremsstrahlung ra- 
diation from electrons in Hll regions, but this becomes dom- 
inant only at high frequencies (> 30 G Hz) and at 1.4 GHz 
only comprises ~ 10% of the radio flux (|Condonlll992l ). 

It is difficult however to explain the linearity and tight- 
ness of the correlation between thermal FIR luminosity and 
non-thermal radio luminosit y. 'Minimum en ergy' estimates 
of magnetic fields in galaxies JBurbidgell 19561 1 H imply a large 
variation between normal galaxies and extreme starbursts 
like Arp220. To explain a constant FIR/radio ratio between 
such disparate systems, complex physical solutions need to 
be provided, for example invoking strong fine-tuning to regu- 
late electron escape and cooling timescales, or short cooling 
timescales with magnetic fields ~ 10 tim es stronger than 
the ' minimum energy' argument suggests ([Thompson et al.l 
I2OO6I). 

IVoelkl (| 19891 ) first suggested a 'calorimeter' model 
whereby both UV light from massive stars and CR electrons 
from SNRs are (a) proportional to the supernova rate, and 
(b) efficiently absorbed and reprocessed within the galaxy, 
so that the respective energy outputs in FIR re-radiation 
and radio synchrotron would both be tied to the supernova 
rate. This theory requires a correlation between the average 
energy density of the radiat i on fie ld and the galaxy magnetic 
field energy density. IVoelkl (|l989h argues this is plausible if 
the origin of the magnetic field is a turbulent dynamo ef- 
fect, since the turbulence would be largely caused by the 



^ Inferring magnetic field strengths of galaxies is problematic 
and depends on many unknowns. The minimum energy argu- 
ment circumvents some of these problems by specifying the mini- 
mum of the total energy density (a function of field strength), 
which occurs at the point where the energy density in parti- 
cles i s in approximate eq uipartition with that in the magnetic 
field iJLongair. M. 3.111994^ . Although there is no physical require- 
ment for this equipartition to occur, it does provide an order-of- 
magnitude estimate for the magnetic field energy density, and 
is physically motivated if one considers that the magnetic field 
is tangled by turbulent motions in the interstellar plasma, and 
that these motions e fficie ntly accelerate cosmic r ay particles (see 
iLongair. M. S.lll994l and iThompson et al]|2006l ). Equipartition 
has been observationally confirmed locally in the Milky Way 
([strong, Moskalcnko & Rcimcr 200C1) and in normal spiral galax- 
ies ([Vallcc 1995; Beck 2000), altho ugh its validity in starbursts 
and ULIRGs is less certain (IThompson et al.ll2006l ). 
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activity of massive stars, and hence correlated with the su- 
pernova ra±e;_AJternativejion-calorirnetric models include 
those of iHelou fc Bicavl () 19931 ). using a correlation between 
disk scale height and the escap e scale length for CR elec- 
trons; and lNiklas fc Beckl |l993), in which the FRC is driven 
by correlations with the overall gas density and equipartition 
of magnetic field and CR energy. 

iBelll ((2003) argues for a 'conspiracy' to diminish both 
the FIR and radio emission originating from star formation 
in low luminosity galaxies when compared with luminous 
~ L* galaxies - without this the relationship would not 
remain linear over th e full luminosity range. Simila rl y, the 
calorimeter model of iLacki. Thompson fc QuataertI ( 2010|) 
invokes conspiracies in low and high gas surface density 
regions to maintain the relationship. The physical origin 
of the FRC therefore is still an open question. A full re- 
view of the theories is beyond the scope of this study, but 
a more detailed discussion of the literature can be found 
in IVlahakis. Eales fc Dunnd (|2007l ). and a m ore in-depth 
treatment is provided in the numerical work of lLacki et al.l 
l|20ld ). 

An ongoing strand of research at the current time is 
the investigation of the FRC at high redshifts and low 
fluxes, in particular whether there is any evolution (e.g . 
Appleton et al..l [20041: iFraver et al.j|2006l: llbar et al .Il200g : 



Garn et al.ll2009l: ISeymour et al.ll2009l : llvison et al..ll2010~ ; 
Sargent et al.j|2010al la). Measurements of any evolution (or 



lack thereof) would improve the accuracy of FIR-/radio- 
estimated SFRs at high redshift, and could shed light on 
the mechanism governing the FRC, as well as highlighting 
differences in the physical and chemical pro perties of star- 
form ing galaxies at high and low redshift (|Sevmour et al.l 
I2OO9I ). 

In the current work we investigate the FRC over a range 
of redshifts, for a sample that is not limited by FIR or radio 
flux. Using Spitzer FIR data and radio data from the Very 
Large Array (VLA) and Giant Metre- Wave Radio Telescope 
(GMRT) , we quantify the FIR-Radio Correlation as a func- 
tion of redshift in massive galaxies selected from a near- 
infrared (NIR) survey of the Extended Chandra Deep Field 
South (ECDFS). We use Equation[T]to define the 'g' index, 
which quantifies the FRC as the logarithmic ratio between 
a monochromatic FIR flux (S^jr, e.g. at 24, 70 or 160^m), 
and 1.4-GHz radio flux (Si/, 1.4 ghz). 



giR = logi 



S. 



Sv,1.4 GHz 



(1) 



We also investigate the effects of using different FIR bands to 
quantify the FRC, and the effects of assumptions about the 
SEDs of the galaxies in the sample. We employ a 'stacking' 
methodology to recover sufficient signal-to-noise ratios on 
faint objects to obtain measurements of the average prop- 
erties of the sample. The data are described in Section [2l 
while the binning and stacking methodologies are described 
in Section [S] The analysis of SEDs and application of K- 
corrections is covered in Section |4l and the results are anal- 
ysed and discussed in Section [S] A concordance cosmology 
of Qm = 0.27, Qa = 0.73, Ho = 71 kms~^Mpc"^ is assumed 
throughout. 



2 DATA 

The ECDFS is a ~ 0.25 deg^ square centred at 
ghggmgQS^ -27°48'20" (J2000). It is a much-studied region 
of sky, with a rich body of published data and studies of ex- 
tragalactic sources at a broad range of wavelengths stretch- 
ing from X-ray to radio regimes. We use radio imaging data 
at 1.4 GHz from the VLA, with a typical rms sensitivity 
across the map of 8/xJy beam"'^ a nd beam dimensions of 
2.8 X 1.8 arcsec (|Miller et al.ll2008h . Imaging data were ob- 
tained at 610 MHz from the GMRT, reaching a sensitivity of 
~_40uJy beam~^ with a 6.5 x 5.4 arcsec beam (|lvison et al.l 
l2010al '). For the FIR, Spitzer MIPS images at 24, 70 and 
160/im were obtained from the FIDEL survey (DR3; Dick- 
inson et al., in preparation). 

To look at a range of galaxy types over a range of 
redshifts we must give careful thought to how the galax- 
ies are selected. For example, selecting radio-bright galaxies 
will naturally favour active radio galaxies, while selection 
at 24/im is likely to favour galaxies with dusty starbursts 
and/or obscured AGN components. These biases will affec t 
the measured value of gm, (see e.g. ISargent et al.l l2010al ). 
There is however a good body of evidence that distant mas- 
sive galaxies in a range of phases of star-formation and nu- 
clear activity can be effectively selected in NIR filters at 
~ 2/im. This part of the spectrum is minimally affected by 
dust absorption, AGN and other components, and hence is 
relatively insensitive to the 'type' of galaxy or the shape of 
its SED. Furthermore it is insensitive to the age of the stel- 
lar population (hence SFR), because the light is dominated 
by old main-sequence stars that make up the bulk of the 
stellar mass in all galaxies. Thus N IR luminosity is primar- 
ily dependent on stellar mass only (jGlazebrook et al.lll995l : 
lGardneillT995[ ). 

In the ECDFS there exist s NIR data from t he Spitzer 
Infrared Array Camera (IRAC: lFazio et aLll2004h . of which 
the two shortest wavelength bands (3.6 and 4.5/im) can 
also be used as an effective tracer of stellar mass (e.g. 
ISerieant et al..ll2008l V A catalogue of IRAC sources matched 
with optical-NIR photometry in the Multiwaveleng th Sur- 
vey by Yale-Chile (MUSYC; ICawiser et aD l20od ) is col- 
lected in Spitzer's IRAC and MUSYC PubUc Legacy of the 
ECDFS (SIMPLE; Damen et al., in preparation). The cata- 
logue was extracted from IRAC 3.6 and 4.5-/^m images, and 
sources for which the mean of the 3.6 and 4.5-/im AB mag- 
nitudes ([3.6] -I- [4.5])/2 < 21.2 were selected, giving a cat- 
alogue of 3841 sources. The source extraction and selection 
is described bv lDamen et al.1 (|2009l ). In SIMPLE the IRAC 
sources have been matched to multiwavelength counterparts 
in the MUSYC catalogue, which contains photometry in UB- 
VRIz'JHK bands. Stars have been identified and excluded 
from the catalogue using the colour criterion J — K < 0.04, 
and potential AGN were rem oved by excluding any matches 
with Chandra X-ray sources (jVirani et al.ll2006l l. Photomet- 
ric redshifts were collated for al l objects in the sample from 
COM BO-17 JWolf et al.j|2004l). and by usi ng the EAZY 
code (Brammcr. van Do kkum fc Coppil 120081 ) as described 
bv lDa nicn ct al. (2009). Damen et al. compared the photo- 
metric redshifts to spectroscopic ones where available, and 
showed that the median {zspec — Zphot) / (1 + z) = 0.033 (0.079 
at 2: Jj 1). As described in Section [3.11 we divide the sam- 
ple into bins with sizes Az/{1 + 2) ~ 0.2 — 0.4, so it is safe 
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to neglect these photometric uncertainties. The final cata- 
logue used in this work contains 3529 sources with photo- 
metric redshifts up to 2; = 2, in the region of the ECDFS 
defined by the rectangle 52°51'48" < RA < 53° 25' 14", 
-28°03'27" < Dec < - 27°33'22" (J2000). 



3 STACKING METHODOLOGY 

3.1 Sample Selection and Redshift Binning 

Since sources are selected by their NIR flux across a range 
of SED types, many are likely to be faint or undetectable 
at the wavelengths of interest, [j In order to probe the evo- 
lution of fluxes as a function of redshift, we stack galax- 
ies into seven bins in redshift and measure median fluxes 
in each bin. The great advantage of this technique is the 
gain in signal-to-noise ratio, as combining many sources 
reduces the random noise while maintaining the average 
level of the signal. This gain is at the expense of knowl- 
edge of the individual galaxies, but with careful appli- 
cation of criteria when binning the galaxies, and with a 
large enough sample, it can reveal properties of galaxies 
below the noise and confusion levels. The technique has 
been used to g reat effect many tira es in the lit e rature ; 
for example b y ISerieant et al..l (l2004h: iDole et al.l JJOO^: 



Bin 



Boundaries 



Median {z) Count N 



l y i&erieant et ai..i iizuu4i i : luoie et ai.i iizuuai: 
ll2007h: iTakagi et al.J (|2007|) : IWhite et all 



llvison et al.,. ___, 

l|200'it l: IPapovich et al.l lJ2007l ): lDunne et al..l (J2009l ). 



We do not know the distribution of fluxes in the stacks, 
but since we select massive galaxies with unknown SEDs 
at a range of redshifts, we may expect to be prone to 
some outliers. For example, radio-bright AGN have unusu- 
ally high radio fluxes and are outliers on the FRC. We can- 
not be certain that these have been successfully removed 
from the sample by cross-matching with the Chandra cat- 
alogue, as we know that there is limit ed overlap between 
X-ray and radio-selected AGN samples ([Rovilos et al.l 120071 : 
IPierce et aLll2010l : iGriflith fc Sternll2010l . although see Sec- 



tion I5.1|l . We therefore used the median statistic to rep- 
resent the properties of the typical galaxies in each stack, 
b eca use unlike the mean, the median i s resistant to out- 
liers JGott et al.ll200ll : IWhite et al]|2007l : ICarilli et al1l2008l : 
iDunneet al.ll2009l '). 

To study redshift evolution, we divided the sample ac- 
cording to the photometric redshifts in the catalogue, and 
stacked each bin into the radio and FIR images. The sample 
was split into the redshift bins given in Table [l] 

The most significant sampling bias that we expect to 
see is that of stellar mass. The stellar masses of gal a xies i n 
the catalogue have be en estimated bylDamen et al.l (|2003 ). 
by SED- fitting with a iKroupal (|2001 ) initial mass function 
(IMF). The effect of Malmquist bias is that the median 
mass in the sample is lower at lower redshifts (because low 
mass galaxies dominate the population), but is higher at 
higher redshifts where only the more massive galaxies are 



^ As an indication of the need to stack, we roughly estimated the 
minimum SFRs detectable at high redshift, using the noise levels 
in Table [2] and calculating SFRs that would be derived from 5it 
detections using the methodology described in Sections 4 and 5. 
These are, from the MIPS and radio fluxes respectively, 870 &i 
13OM0yr-l at 2 = 1, and 4900 & 870MQyr-i at 2 = 2. 
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Figure 1. Scatter plot of stellar masses in the catalogue as a 
function of photometric redshift, with the divisions between the 
redshift bins marked as dashed lines. Large circles mark the me- 
dian mass and redshift in each bin. 



detectable. This is illustrated in Fig. [T] and we test the ef- 
fect on our results in Section (5) 



3.2 Stacking into the VLA and GMRT Radio 

Images 

Radio images have pixel units of Janskys (Jy) per beam, and 
adhere to the convention whereby each pixel value is equal 
to the flux density of a point source located at that posi- 
tion. The assumption that the pixel value at the position 
of each catalogue object gives the correct radio flux density 
for that source is generally good, though requires a small 
correction to give the total integrated flux of the source. 
This integrated-flux correction accounts for any sources be- 
ing extended over more than one beam, and also for any 
astrometric offset between the catalogue coordinates and 
the radio source position. Since we are stacking, it is suit- 
able to consider the overall effect on the median stacked 
source, so the integrated-flux correction is calculated from 
the stacked 'postage-stamp' image of the full sample. This 
image is created by cutting out a 41 pixel squarqj centred 
on each source, and stacking the images by taking the me- 
dian value of each pixel. The integrated flux is calculated 

^ 41 pixels corresponds to 49.2, 164, 328, 20.5 and 61.5 arcscc in 
the 24-/im, 70-/im, 160-/im, 1.4-GHz and 610-MHz bands respec- 
tively. 
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Figure 2. Postage-stamp images of tlie stacked targets in (from left to right) 24/im, 70/^m, 160/im, 1.4 GHz, 610 MHz. Tlie top row is 
the stacked full sample; other rows are the stacks of redshift bins ZB0-ZB6 (from top to bottom) as described in Section l3.1l Images are 
coloured by flux in Jy pixel"^ units (MIPS) and Jy beam~'^ units (radio), and colour scales are included below each. Contours of [3, 5, 
... 17, 19] times the pixel noise arc overlaid. Each tile is a square 41 pixels on a side, which corresponds to 49.2, 164, 328, 20.5 and 61.5 
arcsec for the five bands respectively. 
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using the AlPS package jmfit and the correction is simply 
the ratio of this to the value of the peak (central) pixel in 
the image. 

A higher integrated/peak ratio might be expected in 
the lowest-redshift bins if a large number of resolved sources 
were included in the stack. We therefore tested whether 
the correction varied significantly between different redshift 
bins, by constraining the beam centroid position and ori- 
entation and measuring integrated-/peak-flux ratios in each 
postage stamp. Without these constraints the fitted fluxes 
would be more prone to flux boosting by noise peaks, as 
the signal-to-noise in these postage stamps is low. In both 
1.4 GHz and 610 MHz however, the result for each of the bins 
was consistent (within la) with the result for the stack of all 
sources, indicating that the stacked sources were unresolved 
in all bins. Hence for all bins we used the integrated-flux cor- 
rections from the full stack (which have the smallest errors), 
i.e. 1.55 ± 0.08 for 1.4 GHz and 1.05 ± 0.12 for 610 MHz. 

A further correction would need to be made to point- 
source fluxes measured in the 1.4 GHz image to account for 
bandwidth smearing (BWS), an instrumental effect caused 
by the finite bandwidth of the receiver resulting in sources 
appearing more extended with increasing angular distance 
from the centre of the pointing. Since integrated flux is con- 
served in BWS, the effect is corrected in our data by the 
integrated-flux correction (this explains the large correction 
at 1.4 GHz which would otherwise appear to be inconsistent 
with unresolved sources). 

3.3 Stacking into the Spitzer FIDEL Images 

Measuring fluxes in the MIPS images requires a different 
technique, due to the large point-spread function (PSF) 
which results from the diffraction-limited resolution of 
MIPS. The centre of the PSF can be described by a 
roughly Gaussian proflle, with full widths at half maximum 
(FWHM) of 6, 18 and 40 arcsec in t he 24, 70 and 160-^m 
bands respectively (JRieke et al.ll2004l ). The outer portion of 
the 24- ^m PSF is less predictable and it is known to vary be- 
tween different pointings and different source colourQ For 
this reason, and to allow for potentially resolved sources, we 
chose to measure 24-/xm fluxes by aperture photometry, and 
adopted an aperture of 13 arcsec with an aperture correc- 
tion of 1.16 b ased on the theoret ical STinyTim PReIj and 
the results of iFadda et al..l (|2006h rl to measure total fluxes 
in Jy. Due to poorer resolution in the 70 and 160-/im maps, 
it is sufficiently accurate to measure point-source fiuxes by 
applying a correction to the central pixel value: the factors 
used were 43.04 and 46.86 at 70 and 160^m respectively. 
This converts the fluxes to units of Jy beam~^, and accounts 
for large-scale emission in the wings of the PSF, as well as a 
colour correction. No further correction is required to mea- 
sure total (integrated) fiuxes as we can confidently assume 



at 



'' See MIPS Data Handbook available 

|http : //ssc ■ spitzer ■ caltech ■ edu/mips/dh/ | 
" TinyTim for Spitzer developed by John Krist for the Spitzer 
Science Center. The Center is managed by the California Institute 
of Technology under a contract with NASA. Web page available at 
Ihttp : //ssc . spitzer . caltech . edu/archanaly/contrlbuted/stinyt 
" Sample xFLS and STinyTim PRFs arc available on the SSC 
website: |http : //ssc . spitzer . caltech ■ edu/mips/psf ■ html | 



that none of the sources is larger than the beam in these 
two bands. 

In stacking the FIDEL images it was necessary to ex- 
clude objects close to the edges of the map where the noise 
was higher, to ensure that noise in the stacks reduced as 
l/x/iV and to prevent gradients being introduced into the 
postage stamps. This was achieved by placing lower limits 
on integration time. Limits were chosen based on stacks of 
random positions, resulting in the exclusion of 3.5, 8.6 and 
9.9 percent of the 3529 catalogue sources in the 24, 70 and 
160/im bands respectively. Because these cuts are based on 
integration time alone, there is no correlation with the na- 
ture of the sources themselves, so no systematic effect on 
the measured properties of the galaxies will be introduced. 

Postage-stamp images of the stacked bins in the MIPS 
and radio maps are shown in Fig.[2j including noise contours 
as described in the following section. 



3.4 Analysis of Random Errors in Stacked Flux 
Measurements 

Random errors in any flux measurement arise from noise in 
the image. Simplistically these errors might be expected to 
arise purely from the variance of pixel values in the map, a^, 
and the error on A'^ stacked measurements is then given by 
a l^fN . This assumption is valid for the radio images, so it 
is sufficient to use the rms values at the corresponding posi- 
tions on the rms map. In the MIPS data maps however, pixel 
covariance provides a non-negligible contribution to the er- 
ror, so the rms maps are not sufficient. In order to measure 
the total random error on a measured flux we chose random 
positions in the sample region of the map and selected those 
that fell on empty regions of sky. This was tested by tak- 
ing an aperture of radius 13 arcsec (the radius for aperture 
photometry at 24^m) around each position and measuring 
the standard deviation of pixel values in that aperture. If 
the aperture contained any pixels that deviated from the 
aperture mean by more than 2.5o", then the position was 
discarded. The positions were also required to be separated 
and not overlapping. Thus positions were chosen to repre- 
sent regions of empty sky with no sources. The number of 
positions used was chosen to be 500, to roughly match the 
sample size of bins used for stacking sources. For the 24/im 
case, where aperture photometry was used, the fluxes at the 
500 positions were measured in exactly the same aperture as 
was used for source photometry, and the standard deviation 
of these sky fiuxes was taken to represent the random er- 
ror on an individual aperture measurement. Repeating the 
analysis with different random catalogues of varying sizes 
produced consistent results, and the distribution of fiux val- 
ues in the random catalogues was verified to be Gaussian 
with high certainty. 

The analysis was repeated for the other two MIPS maps 
and also both radio maps, using the central (brightest) pixel 
for flux measurement since this was the method used for 
source photometry in those bands. The resulting error or 
noise values are given in Table [51 The radio error values 
were close to the average value in the rms map, confirming 
that the radio error is equal to the rms value. These values 
im/ iySe'!^}45Miifi'''iltiplied by the aperture correction or integrated- 
flux correction in each case to represent the error on a fiux 
measurement in Janskys. We checked that the noise in a 
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Figure 4. Autocorrelation function W{9) of tlie SIMPLE cata- 
logue. The dashed line is a power-law fit with index -1.15. Error 
bars are twice the Poisson errors, as described in the text. 



Figure 3. Results of measuring noise in stacks of empty sky posi- 
tions. The noise was calculated by stacking TV random positions, 
repeating 500 times, and taking the standard deviation of the 
500 stacked fluxes (see text for details). Results of these stacks 
are shown as symbols, and are in the same units as in Table [2] 
Lines show the least-squares fit to the results of each band with 
a power-law index of —0.5, representing the expected reduction 
of noise as l/y/N for each band. The agreement between results 
and expectations is good. 



stack reduces as roughly l/\/N, as shown in Fig. [S] so the 
errors given in Table [2] are divided by the square root of the 
number of objects stacked, A^. 

In using the median to represent the fluxes of A'^ sources 
in some bin, we must also consider the width of the distribu- 
tion of fluxes in that bin: if this is larger than the estimated 
measurement error/viV then the latter is a poor indicator of 
uncertainty on the quoted median. For this reason we esti- 
m ated la unc e rtaint ies on the median following the method 
of lGott et al.l (|2001f ) and compared them to the estimated 
measurement error in each stack (the value in Table [2] di- 
vided by ViV). At 160^m, 610 MHz and 1.4 GHz we found 
the two to be about equal (see Table [S]), confirming that 
the flux errors we have estimated cover the distribution of 
fluxes in the bins. At 24/im the uncertainty on the median 
was around three times the size of the estimated flux er- 
ror, and at 70^m around twice the size, indicating that in 
these bands the flux distribution in each bin was somewhat 
broader than the estimated errors allowed for. It is possible 
that the method described in preceding paragraphs system- 
atically underestimates the noise in these images, as a result 
of the constraints used to identify empty 'sky' apertures. 
Those constraints were designed to distinguish true read- 
noise on the detector from confusion noise in the sky, but 
the 24-/im image in particular is highly confused, meaning 
that the constraints coud lead to correlation in the empty 
apertures stacked, and increase the chance of underestimat- 
ing the noise. For the analysis of stacked results in this paper 
we therefore quote th e uncertainties on the median following 
the lGott et~al] (|200ll ) method. 



3.5 Background Subtraction and Clustering 

Analysis 

A similar methodology to the random error analysis was 
used to measure the background value to be subtracted. The 
method described above chooses empty apertures containing 
just sky, which is simplistically what needs to be subtracted 
before performing aperture photometry, but in the case of 
MIPS the combination of high source density and low res- 
olution require that source confusion is also accounted for. 
When stacking, the random boosting of fluxes on individ- 
ual sources will average out to a constant correction that 
can be included in the background. Hence when measuring 
the background for subtraction, a catalogue of random po- 
sitions were chosen and stacked, without any criteria on the 
existence or otherwise of sources close to these random po- 
sitions. On average the random catalogue should coincide 
with sources with the same probability as the source cat- 
alogue does. In this context we use the term 'sources' in 
the general sense, meaning any source of flux in the images 
that we stack into, i.e. any object that could boost a mea- 
sured flux at a given position. At this stage we are making 
the assumption that source clustering does not play a part. 
Stacks of 3500 random positions (to match the sample size 
of the source catalogue) were made and repeated 1000 times 
in each of the three MIPS and two radio maps. The mean 
of the 1000 stacked fluxes was taken to be the background 
value, and the standard error was taken to be the uncer- 
tainty; results are given in Table [21 

Any clustering of the sources in the catalogue would 
lead to an increased probability of confusion for a catalogue 
source compared with a random position, hence with in- 
creased clustering the background subtraction becomes in- 
creasingly less effective. In order to estimate the size of this 
effect we would ideally need to understand the correlation 
function of sources in each image on scales smaller than the 
beam size. Since this is not possible, we made the assump- 
tion that the correlation of sources in the images that we 
stack into is approximately the same as that in the source 
catalogue. This may not fully account for confusion if the 
sources in the image are more clustered than the IRAC (cat- 
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Band 


Pixel scale 


PSF/Beam FWHM 


Noise level 


Background level 




arcsec 


arcsec 


MJy 


P'Jy 


24^m 


1.2 


5.9 


62 


-37.60 ±0.04 


70/im 


4.0 


18 


1,200 


+2.2 ±0.8 


160/^m 


8.0 


40 


20,000 


+2, 000 ±10 


1.4 GHz 


0.5 


2.8 X 1.5 


8.83 


-0.014 ±0.005 


610 MHz 


1.5 


7.7 


71.9 


-0.01 ±0.03 



Table 2. Information on the FIDEL and radio images. Measured total noise values represent the la error on a single flux measurement; 
for 24/im this is the noise on a corrected aperture flux, for 70 and 160/im it is the noise on a corrected point source flux, and for the 
radio it is the noise in a beam and does not include the integrated-fiux correction. Background levels are in the same units, these are the 
values subtracted from the median source fluxes. Errors on background fluxes are standard errors from 1000 measurements as described 
in Section 13.51 



alogue) sources, but it does at least remove the possibility 
of double-counting the fluxes of confused sources in the cat- 
alogue. 

We calculated the autocorrelation function W{6) of po- 
sitions in the catalogue, to estimate the excess probability 
of a background source appearing at a radius 9 from a tar- 
get source, compared w ith a random position. We used the 
iLandv fc Szalavl (|l993l ) estimator (Equation[2]) which counts 
pairs within and between the data (D) and random {R) po- 
sitions as a function of annular radius 6. 



WD,D{e) = 



DD-2DR-RR 
RR 



(2) 



Results are shown in Fig. |4j which includes a fit by 
linear regression given by W{9) = 0.000269 6^^'^^ , where 9 is 
in degrees. By dividing the region into four eqaul quadrants 
and comparing the scatter between results in each, we found 
that the standard error was a factor 2.0 larger than the 
simplistic Poisson errors. We therefore quote error bars on 
all correlation functions of twice the Poisson error. 

The strong clustering implies a significant correction to 
the measured fluxes from stacking catalogue sources. The 
correction accounts for the flux contribution from any back- 
ground sources separated by some angular distance 9 from 
the target. This contribution, as a fraction of the average 
source flux, is given by Equation [S] a convolution of the 
correlation function W(9) with the beam profile for the cor- 
responding band (assumed to be Gaussian, exp(— 6'^/2(7^), 
with a =FWHM/(2V21n2)), scaled by the number density 
of background sources n. 



F = n 



W{9)e 



-e^/2o 



2ii9d9 



(3) 



This equation gives the average contribution of con- 
fused sources to a measured flux, hence a correction factor 
of 1/(1 + F) must be applied to stacked fluxes. 

For the 24-/xm case a slightly different convolution must 
be used because aperture photometry is used. The contri- 
bution of a background source to a flux measurement then 
depends not only on where it falls on the beam profile, but 
on how much of its beam falls within the aperture. We com- 
puted the convolution of the 24-^m point-response function 
(PRF jj with the 13-arcsec radius aperture, to give a curve 

'^ We used the empirical xFLS PRF available on the SSC website: 
|http : //ssc. spitzer . caltech ■ edu/mips/psf . html | 



of growth which represents the contribution of a background 
source to the aperture as a function of angular separation 
9. This function is then substituted for the Gaussian beam 
profile in Equation [3] 

This method corrects a stacked flux using the average 
probability of confusion from another source at separation 9, 
scaled by the amount of flux expected from a distance 9 from 
the centre of the beam. When correcting stacks of individual 
redshift bins, we must assume the same level of clustering 
in each bin if we are to use the autocorrelation of the full 
catalogue. To account for the probability of confusion of a 
target from a particular redshift range, while accounting for 
the contribution from background sources at all redshifts, 
we must consider the cross-correlation of the sources in the 
particular range (the 'data' centres, D) with the full cata- 
logue (the 'reference' c e ntres , E\ see Fig.[S|. A modification 
of the iLandv fc Szalavl (|l993r ) method was used to calculate 
the cross-correlation function Wd,e(9), given by Equation|4l 



Wd,e{9) = 



DE-2DR-RR 
RM 



(4) 



The robustness o f the results was tes ted by checking 
against the method of lMasiedi et al.J (|2006l l , which gave in- 
distinguishable results. 

Thus we calculated the average across all the data cen- 
tres, of the excess probability of confusion with any of the 
reference centres. The correction to stacked flux was then 
calculated in the same way as described above, using Wd,e 
in Equation |3] It should be noted that using this estimate 
of the fractional contribution involves the implicit assump- 
tion that the average flux of background sources is equal to 
the average stacked flux. Since we can only correct for con- 
fusion with catalogue sources by this method (i.e. to avoid 
double-counting) this is a reasonable assumption. 

Wc calculated the corrections for three redshift ranges, 
shown in Table O by grouping the bins as follows: 0.0 ^ z < 
0.6; 0.6 ^ z < 1.0; 1.0 < 2 < 2.0. Errors in the table were 
calculated using standard formulae for the propogation of er- 
rors, with the error bars on Wd.e{9) as shown in Fig[S] The 
results in the table are not surprising: 160/im has the low- 
est resolution therefore the greatest confusion, 24/im suffers 
more than 70^m because aperture photometry is used, and 
the radio images have sufficiently high resolution to largely 
avoid confusion. 
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Figure 5. Comparing cross-correlation functions of tliree redshift ranges (D) with the reference catalogue {E) being the full SIMPLE 
catalogue: (a) Wr),E(S) where D is the subset in 0.0 ^ 2: < 0.6; (b) Wd,e(Q) where D is the subset in 0.6 ^ z < 1.0; (c) Wb,e{S) where 
D is the subset in 1.0 ^ z < 2.0. On each of (a)-(c) the black line shows the fit to the autocorrelation function of the full catalogue for 
comparison. Error bars arc twice the Poisson errors, as described in the text. The power-law fits to the four functions are shown in (d), 
where the slope has been fixed by the fit to the full catalogue. 





C 


(cross- correlation) 


C (autocorrelation) 


Band 


0.0< z < 0.6 


0.6<z< 1.0 


1.0 ^ z < 2.0 


All 


24/im 


0.86 ±0.05 


0.74 ±0.04 


0.80 ±0.04 


0.79 ± 0.03 


70/im 


0.90 ±0.04 


0.80 ±0.03 


0.86 ± 0.04 


0.84 ±0.02 


160Aim 


0.81 ±0.05 


0.66 ±0.04 


0.74 ±0.05 


0.72 ± 0.03 


1.4 GHz 


- 


- 


- 


1.000 ±0.001 


610 MHz 


- 


- 


- 


1.000 ±0.002 



Table 3. Correction factors (C) to stacked source fluxes to remove contribution from correlated background sources, using autocorrelation 
of full catalogue and cross-correlations of three redshift ranges with the full catalogue. Errors were calculated as described in the text. 

Corrected flux 5stack,corr = ^stack X C 



4 GALAXY SEDS AND ii'-CORRECTIONS 

4.1 Radio A"— Correction 

Observed fluxes were converted to rest-frame (emitted) 
monochromatic luminosities using Equation [S] which con- 
tains a bolometric iC-correction K{z), accounting for the 
shift of the spectrum in relation to the receiver, and a further 



bandwidth correction [l + z]~^, accounting for the stretching 
of the spectrum in relation to the bandwidth of the receiver 
{(II is the luminosity distance to the source, while z is its 
redshift). 

L^.em = 47rdi S^,ohs K{z) [1 + z]" (5) 

The radio spectrum can be assumed to follow a simple 
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Redshift 

Figure 6. Top: stacked radio fluxes (and la errors) as a func- 
tion of redshift. Note the increase in 610-MHz flux at high red- 
shift which appears to be the main driver for the evolution in 
radio spectral index shown below. Bottom: radio spectral index 
as a function of redshift, calculated from fluxes shown above. The 
dashed line (with error bar) represents the spectral index calcu- 
lated from stacks of the full sample, — 0.74±0.07 (consistent with 
the mean of the indices of all the bins) . The points represent the 
seven redshift bins and are labelled with the corresponding values 
of a, with 1(T errors of ±0.2 calculated as described in Scction [3.4l 
Horizontal error bars mark the widths of the bins. 
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Figure 7. Comparison of four SED templates showing relative 
positions of MIPS filters. Note the differences in the shapes blue- 
ward of the filters, in particular the strength of the PAH and 
silicate features in Arp220, the power-law slope resulting from 
AGN-heated dust in Mrk231, and the cold-dust bump of M51 
peaking at a longer wavelength. 



dian index for each bin to calculate is'-corrections for all 
sources in that bin. 

In Fig. |6]we note an apparent evolution to steeper radio 
slopes at increasing redshift in our sample. A linear least- 
squares fit to the a[z) values gives a slope of — 0.39±0.15; the 
slope is non-zero at the 2.6(t level. This apparent trend is an 
unexpected result, and it is notewort hy that it was not ob - 
served in the stacked 24- /jm sample of llvison et alj (|2010al ). 
who used the same radio data and stacking technique; al- 
though their spectral indices do cover a similar range. The 
possible implications are discussed in Section FS.ll of this pa- 
per, but we note that using a single spectral index of —0.74 
for A" -corrections leads to a slight rise in the q indices in the 
three highest redshift bins (a change of 5q = +0.17 for the 
last bin at (2) = 1.6). 



power law (5*1^ cc v") resulting from the sum of the non- 
thermal synchrotron and thermal bremsstrahlung compo- 
nents; the power law index is ty pically a ~ —0.8 for star- 
forming galaxies (|Condonlll992l ). although stee per indices 
migh t be expected in AGN-dominated sources (|lbar et al.l 
120091 ). The A'-correction to a monochromatic flux with a 
power law spectrum is given by Equation |6l which is inde- 
pendent of the filter transmission function. 



K{z) = [l + z\- 



(6) 



The radio spectral index for each bin was evaluated us- 



ing the stacked fluxes in the two radio bands, 5*^ 
^1.4 GHz jjj Equation [3 (which follows from 5"^ oc v° 
plotted in Fig. [S] 



log [Si 



VSl 



log (610/1400) 



"" and 
and is 



(7) 



These spectral indices were used to iC'-correct each mea- 
sured radio flux using Equation [6l taking the observed me- 



4.2 Infrared Ti"- Correction 

In the mid-/far-infrared part of the spectrum sampled by 
the MIPS bands, the assumption of a simple power law is 
not valid and 7^-corrections must be calculated by evaluat- 
ing Equation |8] which defines the if-correction as the ratio 
of intrinsic luminosity to observed for a general filter trans- 
mission profile T^(i^). 



K{z) = 



!^ n{v)L^{v)dv 



!^T^{u)L^{v[l + z])d{v[l + z]) 



(8) 



This requires knowledge of both the filter transmission func- 
tioijj T^{v) and the SED L„(u). A well-studied local galaxy 
can be used as a template for high-redshift galaxies; com- 
monly used templates in FIR studies include Arp220 and 
M82, which are IR-luminous and therefore considered to be 



* The transmission functions for the MIPS filters are 
available on the Spitzer Science Center website at 
|http : //ssc ■ spitzer ■ caltech.edu/mips/spectral_response . html | 
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K — Corrections to MIPS Luminosities 




Figure 8. A' -corrections for tiie three MIPS bands (from top to 
bottom, 24/im, 70fj,m, 160/^m), calculated from four different SED 
templates. Note the effect of the strong silicate troughs in Arp220 
on the 24- /im flux at redshifts 1.4 (rest-frame wavelength of lO^m) 
and 0.3 (rest-frame 18^m), and the effect of the different dust 
temperatures of the galaxies on the 70 and 160-^m corrections. 



more typical of IR-selected galaxies at high redshift. How- 
ever in the current study there is no reason to presume that 
these dusty IR-luminous objects are representative of the 
stacks, since our sample is selected in the IRAC 3.6 and 
4.5-/im channels. 

We tried a range of different template SEDs to compare 
the results given by the various iiT-corrections. Four tem- 
plates were chosen to represent different types of IR SEDs: 



(i) Arp220 (template from IPope et al..l [ioojl 'l: a bright 
ULIRG with a particularly large mass of hot dust and high 
star formation following a rece nt merger 

(ii) M51 (GRASIL template: [silva et al.|[l99i ) : A typical 
large late-type spiral with moderate star formation and cold 
dust distributed in the spiral arms 

(iii) M82 (GRASIL template: [silva et al.|[l998h : The pro- 



totype hot starburst galaxy, with intense star formation 

probably triggered by a tidal interact ion with M8 1. 

(iv) Mrk231 (GRASIL template: IVega et all l2008h : A 
Seyfert-1, probable merging system, with both starburst 
and AGN components 

The four SED templates are shown in Fig. [T] and the 
_R'-corrections for the three MIPS filters are plotted in Fig.jS] 
The stacking was repeated four times, using a different tem- 
plate for if-corrections each time, to compare the effects 
of different assumptions about the SEDs. The results are 
described in Section [5.21 



5 RESULTS AND DISCUSSION 

5.1 Evolution of Radio Properties of the Sample 

Our results indicate a significant increase in radio luminosity 
with redshift (see Table [Jjl, which is to be expected if the ra- 
dio emission is related to star formation, due to the increase 
in star-formation activity in the most massive galaxies from 
the local universe back to z ~ 2. The apparent evolution in 
radio spectral index over the redshift range (Fig. |6]) is more 
surprising and, notwithstanding the large error bars, hints 
at a fundamental change in the sample demographic, with 
different sources dominating the radio luminosity at z < 1 
and 2 > 1 respectively. The most likely potential contami- 
nant is radio flux from AGN, which would have a different 
spectral index than that from star formation, and would 
also have the effect of boosting radio luminosity. Radio-loud 
AG N source counts a re known to evolve strongly at z > 1 
(e.g. lWalleTaLlbOOSh . 

The effect that AGN contamination would have on 
the median spectral indices is not entirely straightforward. 
While flat {a > —0.5) spectra are ass ociated with radio- 
quiet quasars or low-luminosity AGN (jBondi et al..l 120071 : 
iHuvnh. Jackson fc NorrisI l2007h , steep spectra (a < —1) 
have frequently been used to select powerful radio galaxies 
at high redshift (g enerally z > 2; e.g lde Breuck et al.ll2000l : 
[Pcdani 2003 : Cohcn"eraLll2004l '). This is because AGN radio 
spectra are fiat at low frequencies but steepen at high fre- 
quency, hence steep slopes are observed when the spectrum 
is highly redshifted (the frequency of the turnover varies, de- 
pending on prop erties such as mag netic field strength and 
electron density: iHuvnh et al.l 120071 ). The evolving spectral 
indices seen in Fig. |6]could therefore be a sign of AGN dom- 
inating the radio signal at higher redshifts. 

Matches wit h the Chandra X-ray catalogue 
(jVirani et al.l [200y) have been removed from the sample, 
reducing the likelihood of contamination from unobscured 
AGN. But overla p between X-ray an d radio AGN is known 
to be small (e.g IPierce et al.l |2010| ). In order to identify 
any AGN that are obscured or undetected in X-rays, 
we looked to the MIR fluxes from the IRAC catalogue. 
IRAC colours have been shown to provide some limited 
diagnostics for selecting obscured AGN b a sed on the 
rcst-fram c MIR slope (e.g. [Lacy ct al., 2004; Stern ct ahj 
[2005; Al onso-Herrero et al. J lfoOG: Donley ct al. 2007, 200i). 
In Fig. |9]we plot 58.0Mm/54.5Mm against S2A^im/ Ss.o^im for 
all the objects in the sample. In this plot, the majority of 
objects have 58.0^1x1/54.5^1x1 < 1, well separated from the 
region occupied by AGN at intermediate to high redshifts 
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Figure 9. Obscrvcd-frame Ss.Oiii-n/ S4,.hfj.m - S24tj,ni/Ss.o^ni 
colour-colour plot of all objects in the sample, coloured by red- 
shift bin (for ranges see Table [TJ. Our 24- /^m aperture fluxes are 
used where they are ^ 5cr (roughly 30% of the positions). For the 
positions without a 5cr detection at 24^m, we use the Scr upper 
limit (indicated by arrows). IRAC fluxes are from the SIMPLE 
catalogue. The representative error bar in the top right shows 
the median Icr errors. Some scatter in the vertical axis is intro- 
duced by confusion at 24/im, but positions in the horizontal axis 
are reliable as signal-to-noise is good in both 4.5^m and S.O^m 
and confusion noise is much lower. Tracks of Arp220 and Mrk231 
are overlaid in grey, showing the locus of each on the diagram as 
a function of redshift between and 2. This plot can be used 
as an AGN/starburst diagnostic since AGN have been shown 
to lie to the right, with 5g.o^ni/54.5^ni > 2 at intermediate to 
high redshifts (as the prototype Mrk231 does ) while starbursts 
(such as Arp220) lie above and to th e left (e.g. Ilvison et al.ll2004l : 
IPope et al.]|2008l : ICoDpin et al.ll2010l ^. With the exception of bin 
ZBO (0.0 ^ 2 < 0.4), our data lie well to the left of the plot, in- 
dicating flat spectral slopes at < 8.0/im (compared with Arp220) 
and negligible AGN contamination. The scatter in bin ZBO is 
attributed to strong PAH emission at 7.7^m which makes the 
diagnostic unreliable at low redshift (although it is noted that 
strong PAH emission is generally associated with star formation 
and not AGN). 



(e.g. Ilvison et all |2004| : iLacv et al.ll2004l ). The scatter of 
colours in the lowest-redshift bin is due to the presence of 
the 7.7-fim PAH feature in the sources, and the diagnostic 
is not reliable for this bin. However, the locus of objects in 
all other bins on the diagram clearly indicates that AGN 
contamination is negligible. 

Taking Fig. [9] as evidence against significant AGN con- 
tamination, we deduce that the radio emission originates 
from star-formation activity, and that the index of the non- 
thermal conti nuum from CR ele ctrons evolves. The model 
descri bed by iLacki et al.l (|20ld ) and iLacki fc ThompsonI 
l|201[| ) predicts steep radio spectra as a result of increas- 
ing CR electron losses via inverse-Compton scattering with 
the cosmic microwave background. That model predicts 
these losses to become significant in normal galaxies at 



2 « 2, which seems to be supported by our data. Al- 
ternatively, the steepening spectral slope can be a sign 
of increased electron calorimetry, as electron escape be- 
comes less important relative to electron cooling (including 
synchrotron, bremsstrah lung and inverse-Compton losses; 
iLacki fc ThompsonI [2OI0I ). 

It is possible that the evolution we see in the radio spec- 
tral index is not a variation between galaxies at different 
redshifts, but a function of the rest-frame frequency that 
is observed, and that the assumption of a single power-law 
spectrum is flawed. A curved spectrum with power-law in- 
dex increasing with radio frequency would produce a simi- 
lar effect when vie wed at successively higher redshifts. The 
iLacki et al.l (|2010l ) models predict a steepening of the spec- 
tral index by only ~ 0.05 — 0.1 dex (depending on gas surface 
density) between the frequencies probed by 1.4-GHz obser- 
vations at redshifts from to 2. The evolutionary fit to our 
data indicates a change of —0.8 ± 0.3 over this range, sug- 
gesting that actual evolution with redshift does occur in the 
sample. 

To better understand the change in the sample demo- 
graphic across the redshift bins, we must consider the dis- 
tribution of stellar masses in the respective bins. Malmquist 
bias means that our sample is increasingly dominated by 
the most massive galaxies towards higher redshifts. Stack- 
ing a sub-sample limited to masses > 10^^ M© is unhelpful 
due to small-number statistics, but repeating the stacking 
analysis with a mass limit of log(M*) ^ 10.5 gave results for 
both ct and q that were fully consistent with the full sam- 
ple, although error bars were large (the lowest-redshift bin 
contained too few objects to obtain a reliable value of the 
radio spectral index). Results can be seen in Fig. 1131 



5.2 The Observed and A'— corrected FIR— Radio 
Correlation as a Function of Redshift 

The q index was calculated for each stack and for each FIR 
band i, substituting the monochromatic FIR flux Si (Jy) for 
Sir in Equation [T] Details of the stacked fluxes and q ratios 
can be found in Table[5]at the end of this paper. Fig. llOl dis- 
plays the calculated q values as a function of redshift, both 
before and after if-corrections were applied. Fluxes were 
A'-corrected individually, rather than after stacking, using 
the photometric redshift of each source, and for the radio, 
the stacked spectral index measured in the corresponding 
bin, and for the FIR, the corrections shown in Fig. (8] 

The left panel of Fig. (TDj reveals a very slight down- 
ward trend of q in each of the three MIPS bands, with a 
more bumpy evolution in 524. The anomalies in the observed 
g24 evolution can be explained in terms of the MIR SED, 
which in star-forming galaxies often contains broad PAH 
emission features at 3.3, 6.2, 7.7, 8.6, 1 1.3 and 12.7/jm (e.i 
iRoche et ai]|l99ll : ICenzel et aUll998l : lArmus et all 12007 
These redshifted features can account for the boosting of 
524 in the bins centred at redshifts 0.9 and 1.6. Similarly 
the dip in 524 at redshift 1.3 can be attributed to the broad 
10-/jm silicate trough redshifted into the 24-/xm passband. 
The width of the redshift bins and the use of photometric 
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Figure 10. Far-infrared-radio relation {q) as a function of redshift for the three MIPS bands (from top to bottom, 24/im, 70fj,m, 160^m). 
Ratios of stacked observed fiuxes are plotted on the left. On the right we show ratios of stacked ii'-eorrected fluxes: radio fluxes are 
A"— corrected using the measured spectral index for each bin, and infrared fluxes are X-corrected usin g the four templat es described in 
Section 14.21 Vertical error bars represent estimated Ict uncertainties on the stacked medians, following I Gott et al.l l|200l|): they indicate 
the spread of the data, and are of similar size or larger than the noise, as described in Section 13.41 Horizontal bars in the left-hand panels 
mark the full width of each bin. Horizontal bars are omitted from the right-hand panels for clarity. 



redshifts accounts for the breadth of redshifts over which 
these features appear to have an effect |j 



^ The apparent evolution of q with redshift that is expected from 
various SED templates has been p l otted over observed data by 
sev eral authors inc l uding l lbar et al.l l|2008h . lSevmour et al.l 12009!) 
and lSargent et al.l 1120103^ . In this work we choose instead to plot 



The right panel of Fig. \T0\ shows the effects of K- 
correcting MIPS fluxes using the four SED templates intro- 
duced in Section r4.2l Overall, the M51 template gives rise to 



the X-corrected q- values in an attempt to emphasise the 'excess' 
evolution that may reveal intrinsic changes in the rest-frame flux 
ratios of galaxies within the sample at different redshifts. 
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the least evolution in all three q indices. isT-corrections using 
the Arp220 template give rise to an increasing 524 , due to the 
steeper MIR slope (A < 24/im). In qieo the Arp220, Mrk231 
and M82 templates all exacerbate the downward trend to- 
wards high redshift, while the M51 template removes it, as a 
result of the cooler dust temperature (longer wavelength of 
the peak) in M51. However, none of the templates removes 
the trend in 570, and this could be attributed to a real evo- 
lution in rest-frame flux ratios, or a steeper spectral slope 
at A < 70/im in the galaxies sampled, in comparison to the 
templates chosen. 

5.3 Infrared Spectral Energy Distributions 

Interpretation of our results is evidently subject to the 
assumptions made about the 'average' or typical SED of 
the sample. It is possible to better constrain the FIR K- 
correction by analysing the evolution of MIPS flux ratios 
(colours) as a function of redshift. These colours are sensi- 
tive to the position of the peak of the thermal dust emission, 
hence the temperature of the emitting dust, as well as the 
slope of the SED on the short-wavelength side of the peak. 
Fig. [11] shows the evolution in observed MIPS colours with 
redshift, plotted over the expected tracks for each of the 
SED templates, and reveals that the SED most consistent 
with observed colours at all redshifts is M51. 

The important factor distinguishing the M51 template 
from the others used is the position of the peak of the SED at 
a longer wavelength. M51 is a quiescent star-forming galaxy 
with an IR SED dominated by cold dust, and evidence from 
Fig. [11] therefore points to a cold dust temperature for the 
galaxies in our sample, at least in the first three redshift bins. 
The 70-160-^m colour is directly sensitive to the position of 
the peak at low redshifts, but it is clear from the middle 
panel of Fig. [TT] that over the last four bins the M51 and 
Arp220 templates are barely distinguishable in this colour 
space, so we cannot draw conclusions on the dust temper- 
atures at 2 > 0.8. This is because at these redshifts both 
bands are shortward of the peak of even the hottest IR 
SED, and probe the slope on the Wien side. The 70-160- 
/im colours of the high-redshift bins are consistent with the 
steeper slopes of M51 and Arp220, and not with the shal- 
lower slopes of M82 and Mrk231 (similarly the 24-70-/.im 
and 24-160-/im colours rule out Arp220, due to its strong 
PAH emission). These steeper slopes are potentially an in- 
dication of a stronger contribution from 'cold' dust (in the 
ambient ISM) relative to 'hot' dust (in Hll regions associated 
with star formation) or a dearth of emission from very small 
grains (VSGs); alternatively they could even result from ex- 
tremely optically thick systems where the SED is steepened 
by MIR dust attenuation. In this case, however, we would 
expect to see a stronger 10-^m silicate absorption feature 
such as that evident in the Arp220 iiT-correction at z~ 1.5. 
Our stacked colours are not consistent with such a strong 
absorption which reduces the likelihood that optically thick 
MIR emission is responsible for the steeper rest-frame MIR 
slope at high redshift. 

Cold dust temperatures are nev ertheless consiste n t with 
the conclusions of st udies such as IChapman et al.l (|2005l ) 
and lPope et al.l (l2006l. 2008) for high - redshi ft sub- mm galax- 
ies (S MGs), and ISvmeonidis et all (2003 ). ISevmour et al.l 
(I2OIOI ') and iGiovannoli et al.l (|201oh for 70-/im-selected 




Figure 11. Observed MIPS flux ratios as a function of redsliift 
for the three MIPS bands. From top to bottom: log(S'24/S'7o), 
log(S'7o/5i6o), log(S'24/5i6o)- Overlaid are the tracks of the four 
templates described in Section 14.21 Error bars are as in Fig. 1101 



galaxies at z < 1. There are also parallels with a re- 
cent detail ed study of two ma ssive is'-selected galaxies at 
z ~ 2 by iMuzzin et all (|2010l ). who fitted SEDs to data 
from Spitzer, BLAST and LABOCA jSiringo et al..ll2009l ) 
instruments. Their best fits were star-formation-dominated 
SEDs with I/TiR ~ 1O^^L0, but with cold dust tempera- 
tures, in contrast to ULIRGs in the local universe. Simi- 
lar cool SEDs have also been determined for SMGs out to 
z ~ 1 from BLAST a nd Herschel studies (|Dve et al..ll2009l : 
lAmblard et al.j|2010l '). 

For our sample of massive galaxies, we expect to probe 
the epoch of stellar mass buildup at 2 > 1. Indeed, Tabled] 
shows that both radio- and IR-derived SFRs in our bins 
do reach high values beyond this redshift. It seems a rea- 
sonable assumption that the IR and radio luminosities are 
dominated by star-forming activity, since we do not expect a 
significant contamination from AGN-dominated galaxies in 
the sample (see Section [5. ip . We observe therefore that de- 
spite the tendency towards higher luminosities (and SFRs) 
in the sample at increasing redshifts, there is no evidence for 
a change in the SED towards the templates of local high-SFR 
galaxies such as Arp220 or M82. 

In Table |4] we show indicative TIR luminosities, derived 
for each bin using the rest-frame FIR luminosities in the 
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z Range (2:) 1/1.4 GHz, W Hz 



SFRi 



M 



©yr 



^TIR, L(I 



SFRtir, Mgyr- 



gxiR 



0.00-0.40 


0.21 


(1.03 ±0.13) X 1021 


0.57 ±0.08 


(5.8 ±1.4) X 10^ 


0.61 ±0.17 


2.76 ±0.12 


0.40-0.61 


0.53 


(1.21 ± 0.06) X 10^2 


4.23 ±0.21 


(5.7 ±1.4) X lO^o 


4.8 ±1.3 


2.68 ±0.11 


0.61-0.73 


0.67 


(2.15 ±0.11) X 10^2 


7.50 ±0.37 


(9.9 ±2.5) X lO^o 


8.1 ±2.1 


2.67 ±0.11 


0.73-0.96 


0.87 


(3.15 ±0.19) X 1022 


10.98 ± 0.69 


(1.72 ±0.43) X 10" 


18.4 ±4.7 


2.75 ±0.11 


0.96-1.20 


1.06 


(7.48 ± 0.72) X 1022 


26.1 ±2.5 


(3.27 ±0.82) X 10" 


34.1 ±8.8 


2.65 ±0.12 


1.20-1.42 


1.29 


(1.52 ±0.14) X 1023 


53.0 ±4.9 


(6.2 ±1.6) X 10" 


64 ±16 


2.62 ±0.12 


1.42-2.00 


1.61 


(3.63 ±0.35) X 1023 


127 ±12 


(1.07 ±0.27) X 10^2 


109 ± 28 


2.48 ±0.12 



Ta ble 4 . Stac ked galaxy properties derived from measured fluxes: Rest-frame 1.4-GHz luminosities; SFRs derived from Z/1.4 gHz using 
the iBelll 1120031 ) calibration; Total IR luminosities (Ltir = Lg — lOOOpm ) estimated directly from the if -corrected MIPS fluxes using the 
M51 template; SFRs derived from I/tir using the iBelll 1120031 ) calibration; Corresponding q values calculated as described in the text. 
Errors on -L1.4 gHz are directly translated from the la flux errors. Errors on Ltir are assumed to be 25% as described in the text, while 
errors on SFRs are directly translated from luminosity errors, and do not include any systematics from the conversion to SFR. All SFRs 
are calibrated to a Kroupa IMF. 



MIPS bands, scaled up to I/s-iooo/jm assuming the M51 tem- 
plate. This was done by using the luminosities in the three 
MIPS bands simultaneously to find the best-fitting normal- 
ization of the M51 template. There will be some systematic 
uncertainties in the calibration, and for an idea of the size 
of these we conside r another method to estimate I/tir from 
MIPS luminosities. iDale fc Heloul (|2002l ') offer one such for- 
mula for Lg-iioo/jm, calibrated for a large sample of normal 
star-forming galaxies with a range of morpho logies, colours 
and FIR luminosities (see also lDale et al.ll2001f ): we therefore 
consider it appropriate for M51-like SEDs. Th e uncertainty 
on thi s calibration was shown to be ~ 25% bv lDraine fc Lil 
((20071), and we find that using Dale & Helou's method yields 
values well within 25% (typically 6%, but as much as 16% 
for the highest-redshift bin) of those found using the M51 
template over the same range. Hence assuming systematic 
errors of 25% on Ltir is reasonable. 



Notwithstanding these uncertainties, the results im- 
ply that the typical galaxies sampled have quiescent IR 
SEDs at low z, but rapidly evolve towards higher IR lu- 
minosities at increasing z. By z ~ 2 they appear to reach 
ULIRG luminosities, as star-formation activity becomes sig- 
nificantly more prevalent in m assive galaxies at these red- 
shifts (e.g. iDaddi et al.ll2005bl ). The rise in luminosity with 
redshift which we observe (from ~ 10^*^ to ~ lO^'^Lo) 
may be partially attributed to increasing median stellar 
mass with redshift. This cannot be the full story though, 
since assuming a linear relationship between stellar mass 
and Ltir implies an increase by a factor of 19, whereas 
Ltir increases by a factor of ~ 180, and L1.4 ghz by 
~ 220 over the redshift range. Indeed we know that Ltir 
is linked not to stellar mass itself, but t o SFR, which is 
well-known to rise with increasing redshift (Lilly et al.ll 19961: 



Madau et al.ll996l:[p"erez-Gonzalez et al.l2008l : lDamen et all 
20091 : iMagneUi et al.ll2009l v 



5.4 Evolution in Specific Star Formation Rates 

The S FRs given in Table |4] were ca lculated using t he formu- 
lae oflBelll ((2003 ). which assume a lSalpeteJ (| 19551 ) IMF: 



1.57 X 10~^°Ltir 1 + 



10^ 



SFRtir = <^ 



Ltir 
Ltir > 10" 



1.17 X 10 



-10 r 

ijTIR 



1 + 



10^ 



Ltir 
Ltir sS i6" 



(9) 



SFRi.. 



5.52 X 10" 
5.52 X 10" 



^ r 
-L1.4 GHz, 

2 

rLl.4 GHz, 



L > Lc 
L<Lc 



{ 0.1 + 0.9(L/Lc)0-3- 

(10) 
where SFR is in units of Mqjt ^ , Ltir is given in units of 
L0, L1.4 GHz in W Hz"^ and Lc = 6.4 x 102i W Hz'V These 
conversions were applied to the stacked Ltir and L1.4 ghz 
to obtain SFRs and to stacked Ltir/M* and L1.4 ghz/A^* to 
obtain SSFRs (where M* are t he indiv i dual s tellar masses). 
All SFRs were converted to alKroupal (|200lh IMF by sub- 
tracting 0.2 dex, following iDamen et al.l ( 20091 ) . in order to 
ensure consistency with the stellar masses used. Radio- and 
TIR-derived SFRs appear to be roughly in agreement; the 
TIR values are generally higher although mostly they are 
within the broad error bars given by the calibration of Ltir. 
Agreement naturally depends upon the value of qtir as a 
function of redshift being equal to the local value (e.g. the 
median in Bell's (2003) sample was 2.64). This will be dis- 
cussed in the next section. We note that using a constant 
radio A'-correction based on the overall median spectral in- 
dex of —0.74 reduces radio SFRs in the last three bins but 
does not improve the agreement overall. 

Fig. [12] (black solid points) shows the median specific 
SFRs (SSFRs: calculated source-by-source as the radio 
SFR divided by the stellar mass) as a function of redshift. 
It is immediately apparent that median SSFRs increase 
strongly vnth redshift, indicating a rise in star-formation 
efficiency within the sample at increasing look-back 
times, a res u lt see n many times in the literatur e (e.g . 
Cowie et al.l 19961: iMadau. Pozzetti fc DickinsonI Il998 



Brinchman n fc Elli^j 120 001: iBauer et~ 
2OO5I: [pjr cz- Gonzalez ot all 



12008 



I2OO5I: iFeulner et all 



Damen et all 120091 : lOliver et al.l 



2010bl ) 



Dunne et af] l2009l : 
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Figure 12. Medians specific SFRs from the 1.4 GHz luminosities 
as a function of redshift. SoUd black circles are the medians of the 
full stacked sample; open circles are using the mass limit log M* > 
10.5. Also plotted are the results of lDunne et al.l 1120091 ) with the 
mass limit logM^^ ^ 10.5 (small grey symbols). All vertical error 
bars denote the es timated la uncertainty in the median following 
iGott et al.l ||200U ): systematics in the SFR calibrations are not 
included. Horizontal bars mark the widths of the bins. All SFRs 
in this plot have been converted to a Kroupa IMF as described 
in the text. 



points in Fig. [12] show the results of the mass-limited 
sub-sample (log(Mt) ^ 10.5) in comparison to the full 
sample (filled points). This shows the effect of having lower 
median stellar masses in the full sample at low redshifts 
in particular (because SSFR is a function of stellar mass 
as well as redshift). The consistency between the results 
of this stacking study with the stacked iC'-selected sample 
of lDunne et al.l (|2009l ) seems to support the idea that the 
IRAC selection targets a similar population to K selection. 



5.5 Does the FIR— Radio Correlation Evolve with 
Redshift? 

The M51-corrected q indices in Fig. 1101 which we have 
shown to be the most appropriate, appear to show the least 
evolution in the FRC. It is important to note that the ob- 
served fluxes in the three MIPS bands trace different parts 
of the SED and are affected by different components of 
emission in the source galaxies. In particular, the observed 
160-/im flux between redshifts of - 2 is the closest tracer 
of the FIR peak of the SED (due to emission from large 
graphite and silicate grains, making up the majority of the 
dust mass), and gieo if-corrected with M51 displays no ev- 
idence of evolution. The 924 index is also broadly constant, 
despite 24/im being a closer tr acer of emission from PAHs 
than the bulk of the dust (e.g. iDesert. Boulanger fc Puged 
Il990r ). This lack of evolution is in agreement with previ- 



ous studies of 24-^m-select ed (e.g. Appleton et al.ll2004l ) and 
radio-selected samples (e.g. llbar et al.ll2008l ). Somewhat sur- 
prisingly however, some evolution is still apparent in 570 af- 
ter M51 A'-correction, at around Sa significance. 

Some of the anomalies in the if-corrected q24, and qjo 
graphs could be due to the MIR spectrum and/or the ra- 
dio if-correction. The z ~ 0.9 bin for example coincides 
with the redshifted PAH feature at 12.7^m, and the boost 
in A'-corrected 24-/im fiux at this redshift might be a sign of 
strong PAH emission in the sources. The radio A'-correction 
could also play a part, since in this bin the measured spectral 
index is relatively fiat. This explanation appears likely since 
a similar bump is apparent in gro at the same redshift. Fur- 
thermore, there are particularly low values of 524 and gro in 
the z ~ 1.6 bin, which coincides with the steepest measured 
spectral index. Repeating the stacking analysis using a con- 
stant spectral index of —0.74 ± 0.07 for radio A'-corrections 
was found to have a small effect on both of these bins, chang- 
ing each q index by —0.06 dex at z « 0.9 and -1-0.17 dex at 
z « 1.6. Similarly the values in the intermediate bins at 
z ~ 1.1 and 1.3 were raised by 0.07 and 0.06 respectively 
(changes in the low-redshift bins were negligible); however 
this still leaves a decline of 2cr significance in q-jo when the 
M51 template is used. Clearly it is the FIR SED which dom- 
inates the evolution of monochromatic q indices, and not the 
radio spectrum. 

One factor that could account for this decline in qto 
is a steepening of the continuum slope shortward of 70/im, 
relative to the M51 template. The SED in the MIR region 
(lO/xm < A < 70/im) is thought to be dominated by emis- 
sion from very small grains (VSGs, with radii < lOnm) with 
fiuctuating temperatures resul ting from a mixtur e of ther- 
mal and single-photon heating ([Desert et al.lll99d ) . A steep- 
ening of the the slope shortward of 70/im might be due to 
an increase in the FIR (~ 100/xm) luminosity (dominated 
by big grains) relative to the VSG contribution at shorter 
wavelengths, although this is not clear from the MIPS fiux 
ratios (Fig. fTTj) . 

In t his context it i s inter esting to compare with the 
results of I Seymour et al.l ((20091) , who measured 70-pim fiuxes 
for a sample of faint radio sources and reported a decrease 
in observed gro with redshift (both for detected sources and 
stacks), which is not fully accounted for by the Jf-correction 
of any single model SED. Seymour et al. concluded that 
their stacked data show a discrepancy at 0.5 < 2 < 1.5 
between increasing total Ltir values (estimated from radio 
luminosities) and decreasing 570, implying a change in the 
ULIRG SED at high redshift. Whatever the cause, it seems 
plausible that these two samples are similarly affected. 

Some of the first results from Herschel provide further 
tantalising evidence f or some change in s t ar-for mation activ- 
ity at high redshifts; iRodighiero et al..l (|201Cfl stack ed into 
100 and 160-/xm imaging from t he PACS (jPo glitsc h et al..l 
I2OIOI ) Evolutionary Probe (PEP; iBerta et al.ll2010l '). at the 
positions of IRAC (4.5-/im) sources that were optically clas- 
sified as star-forming and undetected in the 160-/im image, 
divided into bins of stellar mass and redshift. They found 
that SSFRs (derived from IR-I-UV luminosities) followed a 
power-law trend with mass, with an index of —0.25^q\\ 
at z < 1, in agreement with SSFRs fro m radio stacking 
(|Dunne et al.ll2009l ; iPannella et al. J 120091 '), but that the in- 



dex steepened to — 0.50_ 



at 1 < z < 2, deviating from 
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the radio results. A change in the IR SED or gxiR would 
be expected to produce such a deviation between SSFRs de- 
rived respectively from IR and radio (as is suggested by our 
data in Table |4|. 

In Fig. [13] we plot the q indices calculated from L tir 
(listed in Table lU using Eg nation [TT] (|Helou et al.lll985l ): 



gpiR = log 



1/3.75 X 10^ 
W 



_log(i;il^~) (11) 



Here we substitute Ltir for Lfir ~ 1/40- 120^111, (as in iBelll 
I2OO3I and llvison et al.l l2010al . for example) , and this differ- 
ence should be noted when comparing to other work. As an 
indication, the ratio of I/tir/^fir in the M51 template is 
2.1 (which implies gxiR — 5fir ~ 0.32), but this ratio is likely 
to be variable since much of the longer wavelength emission 
can include contributions from dust heated by older stellar 
populations (as discussed for example bv lBcU 2003 ). 

The results for g xm a re sho wn in Fig. [T3] alongside 
the median result of iBelJ (|2003l ) of qtir = 2.64 ± 0.02 
for a FIR-f FUV-selected sample of star-forming galaxies at 
z ~ 0.0. We see that our results are generally a little higher 
than this value at 2 < 1, and it is only due to an apparent 
evolution in our results that they are more in agreement at 
high redshift. The slight discrepancy is just within the er- 
rors allowed by our TIR normalization, and is likely to result 
from a systematic difference in the assumptions made about 
the SEDs and the associated calibration of TIR. 

The slight decline in our qtir values with redshift is 
described by an error-weighted least-squares fit given by 
9TIR oc (1 -|- z)'^ where 7 = —0.11 ± 0.07. Note that stacking 
with the mass limit log(M*) ^ 10.5 gives very similar re- 
sults, fit by the i ndex 7 = —0 . 18 ± 0. 10. In comparison, the 
24/im sample of llvison et al.l (|2010al ) showed evidence for 
evolution over redshifts from to 3, with an error-weighted 
least-squares fit of the same form g iven by 7 = — 0.15±0.03. 
Most recent Iv. llvison et al..l (|2010bl l showed that a sample of 
LIRGS detected by Spitzer and stacked into Herschel imag- 
ing at 100, 160, 250, 350 and 500/im appear to exhibit an 
evolution in qtir over z = — 2, with 7 = —0.04 ± 0.03 (or 
—0.26 ± 0.07, discounting their 16 galaxies at z < 0.5 which 
were poorly matched in Ltir to the higher-redshift bins). 

A slight decline of a similar scale (~ 0.35 dex) in 
Qtir with redshift ( < z < 1.4) was also observed by 
ISargent et al.l (|2010al ). in the median IR/radio ratios of their 
sample jointly selected in the IR and radio. However this 
was at low (2cr) significance and the possibility of intrinsic 
evolution was rejected by the authors because the median 
at 2 ~ 1.4 was within the scatter of their low-z value, and 
moreover because the average at z > 2.5 was very similar 
to the local value. Instead they considered that their sam- 
ple was more contaminated by AGN at increasing redshifts, 
and that the hot dust in these AGN caused 524 ratios to re- 
main constant, while lower abundances of cold dust caused 
570 and qtir to fall. It is interesting to note that we simi- 
larly observe constant 524 and falling 570 and (jttir, but our 
observation of constan t qieo defies a similar explanation. 

In a second paper. ISargent et al.l (|2010a ) extended their 
earlier work using two volume-limited subsets of the joint 
sample: ULIRGs, and sources popul ating the brigh t end o f 
the luminosity function defined by iMagnelli et al.l (|2003 ). 
They showed that for both of these IR-bright populations, 
?TiR was constant out to redshift 1.4. Following a correction 



for increased scatter in their data beyond this redshift, they 
concluded that it remained co nstant out to redshift 2. Qj 
This result disagrees with that of llvison et al.l (|2010al ) , which 
is flux-limited as opposed to volume-limited, showing the 
potential importance of selection effects. 

The decline in our values of qtir could still be caused by 
the same effect that introduces the decline in qjo, since the 
template would underestimate both if the true SEDs were 
steeper at A < 70/im. Alternatively if qtir really declines 
at high redshift then something must be causing galaxies to 
emit less in the IR relative to the radio at increasing red- 
shifts. This could mean either a reduction in optical depth, 
causing more UV photons to escape, or an increase in the 
conflnement and/or reprocessing efficiency of CR electrons 
leading to stronger radio emission. This latter possibility 
cannot be ignored in the light of our observatio n that radio 
spectr al indices steepen at redshifts z > 1, since iLacki et al.l 
(|2010f ) predict that steeper radio spectra are a sign of in- 
creasing electron calorimetry in normal galaxies. 

In spite of these considerations, we remind the reader 
that th e evolution in gTiR is at low significance (similar to 
that of ISargent et al . 2010a), and our data are consistent 
within 1.5(7 with a non-evolution. There is also the poten- 
tial for some bias introduced by the variation of spectral 
index with redshift; applying a constant spectral index of 
—0.74 to the radio iC'-corrections reduces the gxiR evolu- 
tion to a level that is indistinguishable from being con- 
stant: 7 = 0.03 ± 0.07). The evolution in 570 however is 
not fully removed by this change. Using the measured spec- 
tral indices in each redshift bin we fit qro oc {1 + z)'''"' with 
770 = —0.15 ± 0.05, while using the constant spectral index 
we find 770 = —0.10 ± 0.05. Nevertheless, it is emphasised 
that the measured spectral indices should give the most ac- 
curate A'-correction, and Fig. [S] shows that the overall me- 
dian of —0.74 is certainly not appropriate to represent the 
flux ratios in all of the bins. 



6 CONCLUSIONS 

We have studied the FRC as a function of redshift for NIR- 
selected massive galaxies in the ECDFS, a sample which 
is unbiased by star-formation activity. We used a stacking 
analysis to evaluate the ratios of median FIR/radio fluxes 
of all galajcies in the sample, divided into redshift bins. This 
technique traces the typical objects in the population of mas- 
sive galaxies from low redshift back to their formation epoch. 
A thorough analysis of clustering of the sample was used to 
correct for the differential effects of confusion in the three 
FIR bands. iC-corrections were derived in the radio and FIR 
using ratios of observed fluxes, ensuring as much as pos- 
sible a self-consistent analysis. A mass-limited sub-sample 
was also stacked to confirm the robustness of the results to 
Malmquist bias. 

The results for 524, 570 and gieo show a slight decline 
in the observed relations, not dissimilar to the results of 
previous studies, which can be largely accounted for by 



^'^ We note that our data are not affected by the bias in q de- 
scribed by Sargent et al. (2010a, b) due to our selection in the 
IRAC bands. 
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of rec e nt studies inclu di ng Ibar et al. (120081) ; Garn et alj 




2.0 



Figure 13. Total infrared luminosity-to-1.4 GHz luminosity ra- 
tio as a function of redshift for the full sample (black solid points) 
and for the mass limited sample (log(M*) ^ 10.5; grey open 
points). Note that both sets of results decline similarly, imply- 
ing that the evolution seen is not a symptom of Malmquist bias. 
Vertical error bars are dominated by the assumed systematic un- 
certainty of 25% in I/TIR (described in the text), but also include 
the Icr error on the 1.4-GHz flux. Horizontal bars denote the full 
widths of the bins. The star at_z = represents the median 
g.j,jj^ = 2.64 ± 0.02 from iBelll ||2003|1 Our results appear to be 
systematically higher than this value, except where they decline 
at z > 1, potentially as a result of our assumptions about the 
SEDs. 



(l2009l): I Younger et al..l (|2009l ): llvison et all (|2010al ) and 



ISargent et al.l ( 2010al ) . The issue is clouded however by mea- 



surements at 70/^m, which appear to show a declining q in- 
dex with redshift, and when combined into a total IR lumi- 
nosity, likewise show a slight decline (at low significance). 
This most likely implies a steeper spectral slope at wave- 
lengths around 25 — 35/^m (compared with the M51 tem- 
plate), leading to insufficient 70- ^m /('-corrections. But a 
true evolution in the ratios of 70-/im/radio luminosity and 
of TIR/radio luminosity is plausible, considering the appar- 
ent increase in electron-calorimetry behaviour at z > 1, and 
considering the fact that rest-frame 24, 70 and 160-/im fluxes 
can arise from different components o f the dust in a g alaxy . 
It is also co nsistent with the results of lSevmour et al.l (120091 ) 
for gro and llvison et al] (2010a/b) using BLAST/Herschel 
and Spitzer observations to measure gxiR. 

Constraining the FIR SED is one of the greatest prob- 
lems in understanding the FRC and the dust emission in 
general from high redshift star-forming galaxies. Upcoming 
surveys with Herschel, such as the Herschel Multi-t iered Ex- 
tragalactic Survey (HerMES: [oiiver et alTlbOlOal ) and the 
Herschel Astrophysical Te rahertz Large Area Survey (H- 
ATLAS: lEales et al..ll2010l ) are anticipated to revolutionise 
our understanding of these topics by providing deep and 
wide observations spanning the peak of FIR emission across 
the history of cosmic star formation. 



the FIR if-correction using an M51 template. After K~ 
correction qyo is the only monochromatic index to show signs 
of evolution, suggesting that the 70-/im if-correction may 
be less effective as a result of a steep slope in the SED from 
~ 25 — 35/im (corresponding to z ~ 1 — 2) compared with 
M51. 

Observed MIPS colours at all redshifts are more consis- 
tent with the M51 template compared with hotter starburst 
galaxy templates, indicating that the typical IR SEDs of 
stellar-mass-selected galaxies at redshifts up to ~ 0.8 (at 
least) appear to be dominated by cold dust. At higher red- 
shifts it is not possible to constrain the dust temperature 
with MIPS colours, although it is still clear that M51 is 
the closest template. In contrast to this, both radio and to- 
tal IR luminosities rise significantly with increasing redshift, 
as do derived SFRs. Specific SFRs simil arly rise steeply, in 



agreement with results in the literature |Cowic ct al. 



'199^; 



Madau et al.lll998l: [Brinchinann fc EUidfJO OO: Bauer c tld] 
20051: iFeulner et all l2005l: iPerez-Gonzalez et all 12008 



Dunne et al 



Oliver et al.l 



20091: iDamen et al.li2009l: iPannella et al.ll200 



2010bl ). 



The stacked radio data reveal tentative evidence for an 
evolution in radio spectral index across the redshift range, 
an unexpected result that implies some change in the radio 
loss processes in our sample towards higher redshifts. The 
most likely explanation seems to be a shift towards greater 
inverse-Compton losses o f the CR electrons a t z > 1, sup- 
porting the predictions of lLacki fc ThompsonI ((20101). 

Overall our results show evidence that the FRC, mea- 
sured from 24-/im fiuxes or 160-^m fiuxes closer to the FIR 
peak, remains roughly constant up to z ~ 2, corresponding 
to 10 Gyr of cosmic time. This is similar to the conclusions 
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1.2 


21.5 






11.4 


1.2 












610 MHz 


22.2 


3.3 


2.8 


6.8 






20.6 


2.8 






ZBl 


0.53 


528 


24^m 


142.3 


3.1 


9.1 


45.5 


0.94 


0.05 


235.1 


15.1 


1.31 


0.04 








70^m 


1318.9 


53.7 


122.4 


24.5 


1.93 


0.05 


2305.5 


242.4 


2.29 


0.04 








160/^m 


9235.8 


946.8 


1056.6 


9.8 


2.67 


0.07 


4981.4 


674.9 


2.63 


0.05 








1.4 GHz 


14.0 


0.6 


1.0 


23.5 






12.9 


0.8 












610 MHz 


22.8 


3.3 


2.5 


6.9 






19.5 


2.1 






ZB2 


0.67 


529 


24/im 


144.3 


3.1 


10.2 


46.1 


0.97 


0.05 


251.6 


20.1 


1.35 


0.04 








70/im 


1120.2 


53.6 


125.7 


20.9 


1.88 


0.06 


2212.5 


253.3 


2.29 


0.05 








leO/im 


7502.0 


919.8 


694.2 


8.2 


2.60 


0.07 


4152.7 


629.9 


2.59 


0.08 








1.4 GHz 


13.2 


0.6 


0.7 


22.2 






10.9 


0.6 












610 MHz 


20.5 


3.3 


2.8 


6.2 






16.4 


2.2 






ZB3 


0.87 


529 


24/im 


150.5 


3.1 


11.0 


48.1 


1.03 


0.05 


259.2 


21.2 


1.47 


0.04 








70/im 


996.9 


53.4 


107.8 


18.7 


1.87 


0.05 


2262.5 


246.8 


2.37 


0.05 








160fj.m 


6505.4 


932.5 


851.5 


7.0 


2.58 


0.08 


3801.9 


723.8 


2.63 


0.08 








1.4 GHz 


12.1 


0.6 


0.6 


20.3 






9.2 


0.6 












610 MHz 


18.9 


3.3 


4.0 


5.8 






14.2 


2.9 






ZB4 


1.06 


529 


24^m 


126.1 


3.1 


9.0 


40.3 


0.93 


0.06 


262.6 


23.9 


1.32 


0.05 








70/im 


834.9 


53.4 


90.2 


15.6 


1.78 


0.06 


2108.4 


218.9 


2.23 


0.06 








leo^tm 


7843.6 


914.1 


711.8 


8.6 


2.64 


0.07 


5288.1 


625.9 


2.60 


0.07 








1.4 GHz 


12.5 


0.6 


1.1 


21.1 






12.9 


1.1 












610 MHz 


25.7 


3.3 


2.8 


7.8 






25.1 


2.6 






ZB5 


1.29 


265 


24^m 


110.5 


4.4 


9.8 


25.0 


0.77 


0.06 


286.3 


27.7 


1.29 


0.05 








70/im 


741.0 


75.5 


124.8 


9.8 


1.62 


0.08 


1976.6 


354.0 


2.12 


0.08 








160fj,m 


8484.4 


1287.5 


1102.0 


6.6 


2.57 


0.07 


6713.7 


1010.8 


2.64 


0.08 








1.4 GHz 


16.0 


0.8 


1.1 


19.1 






14.7 


1.1 












610 MHz 


33.0 


4.6 


3.7 


7.1 






30.7 


3.4 






ZB6 


1.61 


264 


24/im 


168.2 


4.4 


17.7 


38.0 


0.95 


0.06 


325.3 


24.7 


1.20 


0.05 








70/im 


539.0 


75.9 


101.2 


7.1 


1.48 


0.09 


1489.7 


278.4 


1.81 


0.10 








160/im 


7404.0 


1295.4 


1205.1 


5.7 


2.50 


0.09 


7458.2 


1078.3 


2.56 


0.10 








1.4 GHz 


16.2 


0.8 


1.2 


19.3 






20.3 


1.8 












610 MHz 


40.8 


4.7 


3.5 


8.8 






48.5 


4.2 







Table 5. Summary of stacking results. Columns are as follows: (1) Redshift bin; (2) median redshift; (3) number of objects in stack; (4) 
Band; (5) Med ian observed flux before any corrections; (6) Measured Itr noise (reduced by yN)', (7) Statistical la uncertainty on median 
flux f following iGott et al.ll200lh : (8) Signal-to- noise ratio (Sabs/'^N)] (9) q index for IR band after clustering correction; (10) Error on q 
(using statistical uncertainty from column 7); (11) Median flux following clustering- and X— corrections (using M51 template for MIPS 
and the measured ct{z) for radio); (12) Statistical uncertainty on corrected flux (as in column 7); (13) q index after X-corrections; (14) 
Statistical uncertainty on i^-corrected q. 



